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ABSTRACT 

We combine far-infrared photometry from Herschel (PEP/HERMES) with deep mid-infrared spec- 
troscopy from Spitzer to investigate the nature and the mass assembly history of a sample of 31 Luminous 
and Ultraluminous Infrared Galaxies at z^l and 2 selected in GOODS-S with 24 fjrn fluxes between 0.2 
and 0.5 mJy. We model the data with a self-consistent physical model (GRASIL) which includes a state- 
of-the-art treatment of dust extinction and reprocessing. We find that all of our galaxies appear to require 
massive populations of old (> 1 Gyr) stars and, at the same time, to host a moderate ongoing activity of 
SF (SFR < 100 M0/yr). The bulk of the stars appear to have been formed a few Gyr before the observa- 
tion in essentially all cases. Only five galaxies of the sample require a recent starburst superimposed on 
a quiescent star formation history (SFH). We also find discrepancies between our results and those based 
on optical-only SED fitting for the same objects; by fitting their observed Spectral Energy Distributions 
with our physical model we find higher extinctions (by AAy ~ 0.81 and 1.14) and higher stellar masses 
(by ALog(M*) ~ 0.16 and 0.36 dex) for z~l and z~2 (U)LIRGs, respectively. The stellar mass difference 
is larger for the most dust obscured objects. We also find lower SFRs than those computed from Xir 
using the Kennicutt relation due to the significant contribution to the dust heating by intermediate-age 
stellar populations through 'cirrus' emission (^73% and ~66% of total Ltr for z ~ 1 and z ~ 2 (U)LIRGs, 
respectively). 

Subject headings: galaxies: evolution — galaxies: general — galaxies: interactions — galaxies: starburst 
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1. Introduction 

The spectral energy distribution (SED) of a galaxy 
contains valuable information about its physical prop- 
erties, including the stellar, gas and dust content, the 
age and abundance distribution of the stellar popula- 
tions resulting from the star formation history (SFH), 
and their interaction with the interstellar medium 
(ISM). The study of the SED therefore offers the most 
direct way to investigate galaxy formation and evolu- 
tion, both through direct observations and correspond- 
ing theoretical modelling. The different physical pro- 
cesses occurring in galaxies all leave their imprint on 
the global and detailed shape of the spectrum, each 
dominating at different wavelengths. Therefore, by 
analysing and predicting the whole spectral range one 
can hope to deconvolve and interpret all the informa- 
tion contained in the SED in terms of the SFH and 
galaxy evolution in general. 

Stellar sources mainly emit in the UV/optical to 
NIR spectral range, and the SED in this wavelength 
region is therefore heavily influenced by the star for- 
mation history of the galaxy. The SED from a few 
fjum to the sub-mm (the IR region) is dominated by 
the interaction of dust grains with stellar radiation. 
Basically dust absorbs and scatters photons, mostly 
at wavelength A < 1 fim and thermally emit the ab- 
sorbed energy in the IR. The resulting SED is often 
substantially changed and in many relevant cases rad- 
ically modified by the presence of dust. 

A detailed modelling of the entire SED of galaxies 
is therefore very complex. Several different approaches 
have been proposed, depending on the purpose of the 
applications. Some works (e.g. Devriendt, Guider- 
doni & Sadat 1999; Chary & Elbaz 2001; Dale et 
al. 2001; Dale & Helou 2002; Galliano et al. 2003; La- 
gache, Dole & Puget 2003; Burgarella et al. 2005, 
da Cunha, Chariot & Elbaz 2008) have proposed 
semi-empirical treatments of the SEDs. The aim 
of these works is in general to interpret large sam- 
ples of data, requiring fast computing times making 
use of observationally or physically motivated SEDs. 
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There are then other codes, as Hypcrz (Bolzonella 
et al. 2000) or LePhare (Arnouts et al. 1999) which 
usually perform a template-based SED fitting pro- 
cedure based on a standard % 2 minimization method. 
Usually here dust effects are accounted for by comput- 
ing the attenuation of the light from a source placed 
behind a foreground screen of dust and assuming dif- 
ferent extinction curves. Other works are based, in- 
stead, on theoretical computations in order to have 
a more general applicability in terms of interpreta- 
tive and predictive power. Several papers deal with 
the radiative transfer (RT) in spherical geometries, 
mainly aimed at modelling starburst galaxies (e.g. 
Rowan-Robinson 1980; Rowan-Robinson & Craw- 
ford 1989; Efstathiou, Rowan-Robinson & Siebenmor- 
gen 2000; Popescu et al. 2000; Efstathiou & Rowan- 
Robinson 2003; Takagi, Arimoto & Hanami 2003a; 
Takagi, Vansevicius & Arimoto 2003b; Siebenmorgen 
& Krugel 2007; Rowan-Robinson 2012). Early mod- 
els of this kind did not include the evolution of stellar 
populations. Silva et al. (1998) were the first to couple 
radiative transfer through a dusty ISM and the spec- 
tral (and chemical) evolution of stellar populations. 

Modelling the emission from stars and dust consis- 
tently in order to get reliable estimates for the main 
galaxy physical parameters, like stellar mass (M*), 
star formation rate (SFR) and the average extinction 
(Ay), involves solving the radiative transfer equation 
for idealised but realistic geometrical distributions of 
stars and dust, as well as taking advantage of the full 
SED coverage from UV to sub-mm. The full SED al- 
lows the total luminosity of the galaxy to be robustly 
constrained without relying on extrapolations based 
on optical data alone. While the shape of the optical- 
UV SED indicates some level of dust extinction, there 
may be stellar populations (typically the newly born) 
completely obscured by dust. The only way to over- 
come this problem and quantify how much stellar light 
(and therefore mass) is missing from the optical, is to 
consider the dust-absorbed galaxy luminosity through 
the dust re-emitted spectrum in the IR, and model it 
together with the optical spectrum with a complete 
radiative transfer method. 

With total infrared luminosities between 10 11 — 10 12 
Lq and > 10 12 L©, respectively, Luminous and Ultra- 
luminous Infrared Galaxies [(U)LIRGs hereafter] are 
among the most luminous and complex extragalac- 
tic objects we can conceive, including all varieties of 
young and old stellar populations, dust absorption, 
scattering, grain thermal re-radiation, and AGN emis- 



2 



sion (Lonsdale ct al. 2006). Although they are quite 
rare in the local universe, they dominate the cosmic 
star formation rate and the FIR background at z > 

1. Therefore they are suitable laboratories to study 
the main physical processes which drive galaxy forma- 
tion and evolution. Observations of local ULIRGs have 
shown that they are dominated by strong interactions 
and mergers and the fraction of mergers/interactions 
among them has been found to be strongly correlated 
with their IR luminosity, such that lower luminosity 
LIRGs (Lir < 10 115 L Q ) are ordinary disks while the 
highest luminosity ULIRGs are advanced stage merg- 
ers. 

High-z LIRGs and ULIRGs seem to be, instead, 
not equivalents of ULIRGs in the nearby universe, but 
rather upscaled versions of normal galaxies (Elbaz et 
al. 2011; Nordon ct al. 2012; Symeonidis et al. 2009; 
Rujopakarn ct al. 2011). An implication is that high- 
z ULIRGs could well have gradually evolving SFHs, 
and be not necessarily associated with merger events 
as they are at z^0. 

This paper is the first of a series dedicated to 
a detailed physical investigation of the spectro- 
photomctric properties of dusty high-rcdshift (U)LIRGs 
in an attempt to achieve a deeper understanding of 
these sources. 

Here we concentrate on a small sample of high-z 
(U)LIRGs with the currently richest suite of photo- 
metric and spectroscopic data, available for the first 
time at redshift z ~ 2, combining deep Herschel imag- 
ing with ultra-deep IRS spectra from Spitzer. The 
data is analysed with a state-of-art chemo-spectro- 
photometric model (GRASIL: Silva et al.1998) includ- 
ing self-consistent treatment of dust absorption and 
reprocessing based on a full radiative transfer solution. 

In section Sj2] we describe our data sample. The 
physical modeling performed for this study is pre- 
sented in section $3] The results concerning our SFH 
and SFR estimates are discussed in section i]4] While 
the predicted stellar masses and average extinctions 
for these high-z (U)LIRGs arc presented in section §5\ 
The conclusions can be found in section 

Our cosmological model assumes Ho=71 km 
Mpc-\ S1 A =0.73, H M =0.27. 

2. The (U)LIRG data sample 

Our high-z (U)LIRGs have been selected from the 
sample of 48 IR-luminous galaxies in GOODS-S pre- 
sented by Fadda et al. (2010, F10 hereafter). It in- 



cludes the faintest 24 (im sources observed with IRS 
(S24 ~ 0.15 — 0.45 mJy) in the two redshift bins 
(0.76 — 1.05 and 1.75 — 2.4) and samples the major 
contributors to the cosmic infrared background at the 
most active epochs. At these redshifts, the 24 fjrn 
fluxes translate to infrared luminosities roughly in the 
ranges for LIRGs at z ~ 1 and ULIRGs at z ~ 2. 
This sample is therefore crudely luminosity selected 
(with the S24 limits corresponding to different lumi- 
nosity ranges for the z^ 1 LIRGs and z^ 2 ULIRGs) 
and F10 did not apply any other selections. 

These objects have been selected by F10 mostly 
using the photometric redshifts computed by Caputi 
et al. (2006) since only 50% and 10% of them had 
spcctroscopic-z in the z^ 1 and 2 bins, respectively. 
As emphasized by F10, except for having higher dust 
obscuration, these galaxies do not have extremely de- 
viant properties in the rest-frame UV/optical com- 
pared to galaxies selected at observed optical/near-IR 
band. Their observed optical/near-IR colors are very 
similar to those of extremely red galaxy (ERG) popu- 
lations selected by large area K-band surveys. 

In this paper we concentrate on galaxies powered 
by star formation. Therefore we excluded from our 
present analysis objects previously classified by F10 as 
AGN-dominated on the basis of several indicators such 
as broad and high ionization lines in optical spectra, 
lack of a 1.6 /mi stellar bump in the SED, X-ray bright 
sources, low mid-IR 6.2 /im EW and optical morphol- 
ogy (see also Pozzi et al. 2012). We also excluded 
the sources falling outside the multiwavelength MU- 
SIC catalogue (Santini ct al. 2009), which is our refer- 
ence for complementary optical/near- IR photometry. 
Among the 31 (U) LIRGs fulfilling these requirements, 
10 are at z ~ 1 (all in the nominal LIRG regime) and 
21, mostly ULIRGs, are at z ~ 2. The GOODS-S field 
was observed in 2010 with both Herschel (Pilbrat et al. 
2010) SPIRE (Griffin ct al. 2010) and PACS (Poglitsch 
et al. 2010), giving a total of six bands from 70 to 500 
/im. SPIRE and PACS data are taken, respectively, 
from the Herschel Multi-tiered Extragalactic Survey 
(HerMES: Oliver et al. 2012) and the PACS Evolution- 
ary Probe (PEP: Lutz ct al. 2011) programs. Typical 
noise levels are ~ 1 mJy for PACS 70 — 160 /im and 
~ 6 mJy for SPIRE 250 — 500 /im, including confusion. 

3. Modeling the SEDs of high-z (U)LIRGs 

A physical characterization of the ULIRG phe- 
nomenon requires a multi- wavelength approach and a 
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Par. 


Unit 


Description 


class 


range 




VSch 


Gyr- 1 


SF efficiency 


SFH: Quiescent 


[0.3- 


2.0] 


Tint 


Gyr 


infall timescale 


SFH: Quiescent 


[0.01 - 


-3.0] 


%M b 


— 


% of gas mass involved 


SFH: burst 


[0.1- 


15] 






111 LI1C UUISl ai tb urs t 










Myr 


burst e-folding time 


SFH: burst 


[10 - 


80] 


y*moi 




fraction of gas in MC 


basic 


[0.05 - 


-1.0] 






with respect to the diffuse comp. 








n 




MC optical depth at 1 /jm 


basic 


[5.0- 


108] 


^csc 


Myr 


escape time of newly born 


basic 


[1 - 90] 






stars from their parent MCs 








P 




sub-mm dust spectral slope 


basic 


[1.5- 


2.0] 




kpc 


core radius 


geometrical: King's profile 


[0.01 - 


- 2.0] 



Table 1: Main GRASIL input parameters and their range of values (see Silva ct al. (1998) for a full description). 



detailed treatment of the effects of dust. The spectro- 
photometric + radiative transfer code GRASIL (Silva 
et al. 1998; 2011) satisfies these requirements. It 
computes the SED of galaxies from far-UV to radio, 
with a state-of-the-art treatment of dust reprocessing. 
It includes the radiative transfer effects of different 
dusty environments - star forming molecular clouds 
(MCs), diffuse dust (cirrus), dusty envelopes around 
AGB stars, and a geometrical distribution of stars 
and dust in a bulge-like (King profile) and/or a disk 
(double exponential) profile. The accounting of the 
clumping of young stars and dust within the diffuse 
medium gives rise to the age-dependent attenuation. 
The dust model consists of grains in thermal equilib- 
rium with the radiation field, and small grains and 
PAH molecules fluctuating in temperature. 

The input star formation histories are computed 
with CHE_EVO (Silva 1999), a standard chemical evo- 
lution code that provides the evolution of the SFR, 
M gas and metallicity, assuming an IMF, a SF law 
SFR(t)=z/Sch • M gas (t) k + f(t) (i.e. a Schmidt-type 
SF with efficiency ^sch, and a superimposed analyt- 
ical term to represent transient bursts possibly re- 
lated to a galactic merger, and an exponential infall 
of gas (dMi n f/dt oc exp(— t/r in f)). A very short infall 
timescale, r; n f, can be used to have a so called 'close 
box' chemical evolution model, which ensures that the 
gas going to form the galaxy is all available at the be- 
ginning. This short value of the infall timescale allows 
us to accrete enough mass in stars of different ages to 
reproduce the shape of the UV-opt-NIR SED. In our 



library of SFHs we consider a wide range of values for 
Ti n f from very low (0.01 Gyr) to 2-3 Gyr the latter 
producing longer phases of gas accretion. In the fol- 
lowing we have adopted k=l, f(t) exponential, and the 
Chabrier IMF. Our reference library of SSPs is from 
Bressan et al. (1998, 2002), which directly include the 
effects of dusty envelopes around AGB stars. 

We have used a large and fine grid of theoretical 
SEDs, generated by CHE.EVO+GRASIL. Our SEDs 
span a wide range of input parameters: star forma- 
tion history; obscuration times; dust opacities etc (see 
Tab. [1}. This grid is first used to explore the param- 
eter space. Then an object-by-object analysis of the 
(U)LIRG SEDs was performed in order to refine the 
fit. We reproduced the observed broadband SEDs of 
our 31 (U)LIRGs, including a fit to the IRS spectra, 
which are useful to constrain the physical state of the 
warm ISM and the PAH intensity. 

Figure [1] shows the median best- fit SFHs for the 10 
z ~ 1 LIRGs (top) and 21 z ~ 2 (U)LIRGs (centre and 
bottom). The insets report the mass in living stars and 
the mass of interstellar gas. 

As discussed later, the fitting process appears to 
constrain several of the parameters ruling the SFH, 
in particular T- ln { and JAjch- Small values for r in f and 
high values for vsch , i.e. an early fast and efficient SF 
phase as in the centre panels, were required for 16/31 
objects, 4 LIRGs and 12 z-~2 (U)LIRGs, showing the 
strongest stellar bump in the rest-frame near-IR (as, 
for example, the ones in Fig. [2] top-center and bottom- 
left, right). Smoother SFHs with longer r- m f (Tj n f > 
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Fig. 1. — Median SFHs of z <~ 1 (top) and z ~ 2 (centre and bottom) (U)LIRGs modeled with continuous SF (left) 
or a starburst (right). The insets report the time evolution of the mass in living stars and the mass of interstellar 
gas. The median values for T- m f are 0.5 and 0.28 Gyr in the top-left and top-right panels respectively while they are 
much shorter and equal to 10 Myr for the centre panels. For the nine z~2 (U)LIRGs requiring the SFHs represented 
in the bottom panels the median r- m [ are 1 Gyr for both. 
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Fig. 2. — Typical GRASIL best-fits (solid black line) to the observed SED (red circles) of three z <~ 1 (top) and 
z ~ 2 (bottom) (U)LIRGs. IRS spectra appear in the inset window (magenta line). The color-coded lines represent 
the unextinguished starlight (orange dot-dashed), extinguished starlight (green dashed), cirrus emission (magenta 
dotted) and MC emission (blue long dashed). 
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Fig. 3. — Left: Comparison between the best- fit xt relative to the two different prescriptions for the SFH (continuous 
vs starburst) for our high-z (U)LIRGs. For 4/10 z ~ 1 LIRGs (shown as squares) a starburst SFH gave 'acceptable' 
best-fits. For the remaining 6 LIRGs a starbursting SFH was rejected with high-confidence (Axl > 1). Label 
numbers indicate the object IDs. Right: Comparison of our GRASIL-estimated SFRio with the SFR derived from 
the 8 to 1000 /im luminosity Lir using the Kennicutt (1998) calibration, SFRjt. Filled and open black squares and 
circles represent the different SFRk derived by considering the total Ljfj (cirrus + star forming MCs) and the Lmc 
(MC only), for our z ~ 1 and z ~ 2 (U)LIRGs respectively. Green filled squares and circles are z ~ 1 and z ~ 2 
(U)LIRGs for which a burst is not required, while red filled circles are z ~ 2 ULIRGs requiring the presence of a 
moderate burst from an application of the F-test. Crosses are local (U)LIRGs while magenta stars are high-SFR 
galaxies discussed by Rodighiero ct al. (2011, Rll hereafter). Datapoints with errorbars are median values with 
associated semi-interquartile ranges. 



0.2 Gyr and T m f > 1 Gyr) are instead required for the 
z~l LIRGs and 9 z^2 (U)LIRGs shown in the top and 
bottom panels, respectively. These galaxies present al- 
most 'flat' rest-frame NIR bands and higher UV fluxes. 
Figure [2] summarizes the typical best-fits SEDs of our 
sample showing as an example three z ~1 LIRGs (top) 
and three z -2 (U)LIRGs (bottom). The best-fit SEDs 
of all our galaxies are shown in the Appendix. 

4. Evaluating the Star Formation Rates and 
SFH 

Given the detailed shape of the broadband SED, 
our physical analysis seems to be able not only to give 
an estimate of the instantaneous SFR but also to give 
important hints about the main parameters ruling the 
source's past SFH, i.e. T; n f and vsch as shown in fig. 



ru 

We have investigated here both SF models with and 
without a starburst on top of the Schmidt-type part of 
the SF law (see fig. [1]). We hereafter refer to these as 
the starburst and continuous models, respectively. For 
the majority of our (U)LIRGs, a suitable calibration of 
the Ti„f and vsch allowed us to obtain good fits to the 
observed SEDs with the continuous models. Figure [3] 
(left) compares the best-fit xt relative to these two 
SFHs considered for our high-z (U)LIRGs. 

If we consider first the more extreme cases, i.e., our 
z ~ 2 ULIRGs, for 5 out of 21 sources (red filled circles 
in Fig. [3]) , an application of the F-test to our \ 2 analy- 
sis requires the presence of a moderate starburst. Even 
in these five objects, however, the gas mass involved in 
the SB amounts to only a small fraction of the galactic 
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mass, < 4%, and all of them are observed just at the 
end of the burst event. A peak-phase SF usually pro- 
duces intense far-IR spectra, with MC emission domi- 
nating the MIR and FIR spectrum. Conversely a more 
continuous SFH generates colder far-IR spectra with 
cirrus emission in the sub-mm and larger PAH contri- 
butions around 10 /an. Two examples are shown in 
Fig. [5j This figure shows the typical contributions to 
the composite SED of a highly star forming and grad- 
ually evolving galaxy (left), where the cirrus emission 
dominates the IR/sub-mm wavelengths, and a pure 
starburst (right) whose IR fluxes are dominated by 
MC emission. Both galaxies are taken from the sam- 
ple of high-z PACS-detccted off-main-sequence sources 
(the so called "outliers") of Rll and are modeled with 
GRASIL. This sample includes all the objects (~ 50) 
with sSFR a factor of ~ 8-10 higher than the one of 
MS galaxies. The ~ 12 galaxies shown as magenta 
stars in Fig. [3J represent the most extreme cases with 
the highest SFRs. For the remaining 16 sources of our 
sample at z ~ 2, a continuous model gives a perfect ac- 
count of the data, although without excluding a small 
starburst contribution (black filled circles in Fig. |3|). 

As illustrated in Fig. [TJ in most of our objects the 
only significant effect of introducing a SB event is to 
slightly increase the age at which the galaxy is ob- 
served, i bs, changing from 4.5 to 5 and from 2 to 3 
Gyr from left to right, for z ~ 1 and z ~ 2 (U)LIRG 
respectively. The effect of an ongoing SB on the galaxy 
SED is mimicked, in a continuous SFH, by considering 
the epoch of observation to be closer to the peak phase 
of SF. 

For the z ~ 1 LIRGs, all objects are consistent, on 
the basis of the F-test application, with a continuous 
SFH. 

Table [5] shows the main best-fit physical parame- 
ters derived from our analysis and their mean values 
and standard deviations. We report here three differ- 
ent estimates of SFR, SFRk representing the SFR de- 
rived from the 8 to 1000 fim luminosity Im using the 
Kennicutt (1998, K98 hereafter) calibration, and our 
best-fit SFRs averaged over the last 10 (SFRio) and 
100 (SFRioo) Myr. Fig. [3] (right) shows a comparison 
of our SFRio with SFRk- The latter is defined in the 
limit of complete dust obscuration and dust heating 
fully dominated by young stars (see K98). 

Indeed it assumes that the Lboi of a constant SF 
lasting 100 Myr is totally emitted in the IR (K98; 
Leitherer & Heckman 1995, LH95 hereafter). For a 
constant SF, the Lboi after the first 10 Myr evolves 



relatively slowly because the rate of birth and death 
of the most massive stars (with lifetimes < 10 Myr 
and dominating the Lboi) reaches a steady state (see 
Fig. 2 and 8 of LH95). The K98 SFR/L Fm calibration 
adopts the mean bolometric luminosity for a 10-100 
Myr continuous SF, solar abundance, Salpeter IMF 
(which we have rescaled to the Chabrier one) of the 
starburst synthesis models of LH95, and assumes that 
LFiR=Lboi- Therefore we consider SFRio as the best 
indicator of the current SFR to be compared to this 
calibration. 

As apparent in Fig. [3] (right), our inferred SFRio for 
these high-z (U)LIRGs are systematically lower than 
those based on the K98 calibration (filled circles and 
squares). This is due to the significant contribu- 
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Fig. 4. — Correlation between the SFRio/SFRk ra- 
tio and the fractional contribution of MC emission to 
the total Lir. Points are the same as in fig. [3] Here 
we show the eight local (U)LIRGs for which the frac- 
tional contributions of MC and cirrus emission have 
been already computed by us. Black stars are the two 
representative cases of Rll off-MS star forming galax- 
ies shown in fig. [5] We can see here that on average 
galaxies with higher SFRio/SFRk ratios have total IR 
luminosity mostly contributed by MC emission with all 
of our (U)LIRGs lying on the lower relation. Of course 
the scatter we see in this figure as well as in fig. [3] is 
due to the wide range of ages of stellar populations. 
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Fig. 5. — Figure shows two representative cases of the Rll outliers defined in section|H Left: the galaxy is modeled 
as gradually evolving spheroid and have total Lir dominated by cirrus emission. Right: the galaxy is dominated 
by an ongoing burst of star formation in MCs. While the galaxy on the right has SFR in agreement with the one 
estimated from Lir using the Kennicutt calibration, the one on the left follows our calibration. 



tion of cirrus emission to the total Lir (<~ 73% and 
~ 66% for z~ 1 and z^ 2 (U)LIRGs respectively) 
whose heating source includes already evolved stellar 
populations (ages older than t osc ). As cautioned by 
K98, if all of this energy is ascribed, through this cal- 
ibration, to the recent SF, the SFR is overestimated. 
Instead, given the characteristics of the K98 calibra- 
tion, a fairer comparison between SFRio and SFRk, 
when an ongoing starburst is not clearly dominating, 
is to assign to the latter the Lir by MC only (empty 
cirles and squares in Fig [3] (right)). In this way most 
of the points spread around and above the K98 line. 
The spread has to be ascribed to a range of different 
and non constant SFR and SFH, and t eS c (ranging be- 
tween a minimum of 3 Myr and a maximum of 90 Myr, 
with median values corresponding to 14 and 6 Myr for 
z~l and z^2 (U)LIRGs respectively), i.e. to a spread 
of stellar ages contributing to Lir. Averaging out our 
results for LIRGs and ULIRGs we find the following 
calibration between the total IR luminosity and SFR: 

SFR[M Q /yr] ~ (5.8 ± 0.4) x 10" 11 L m /L & , (1) 

approximately a factor 1.7-2.5 different from the clas- 



sical Kennicutt relation, both scaled to a Chabricr 
IMF. 

In Fig. [3] (right) we also plot the SFR i0 vs SFR K 
(from total Lir) derived with the same physical analy- 
sis, for a sample of local (U)LIRGs (crosses) and some 
of the most highly star forming objects of Rodighiero 
et al. (2011) (magenta stars). 

The galaxies shown in this plot (excluding the open 
points which are the same galaxies as the filled ones, 
but considering only the MC contribution) seem to 
define a broad and continuous distribution contained 
between two correlation lines: an upper envelope dom- 
inated by objects with an ongoing burst of SF in MCs 
as for several local (U)LIRGs and peak-phase SBs, and 
a lower one mostly populated by secularly evolving 
galaxies and late starbursts (as the 5/21 objects of 
our sample discussed above), with substantial contri- 
butions of cirrus emission to the total Lir. 

Indeed the relative contributions of MC emission to 
the total Lir for these sources are shown in Figure|4]as 
a function of the SFRio/SFRr ratio: galaxies falling in 
the lower envelope, with low values of SFRio/SFRr — 
0.5, have small contributions by MC, while those with 
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Table 2: Estimated and average values of the main physical parameters, M+, Ay and SFR, derived from the GRASIL best-fits to the 31 
(U)LIRGs, compared with values based on HYPERZ. 



higher SFRio/SFRk values are dominated by MC 
emission. 

An illustration of this is reported in Fig. [5] for the 
two PACS-sclccted sources of Rll: the one dominated 
by MC corresponds to the starred object in the upper 
figure, while that one dominated by cirrus emission 
falls in the lower part and shows minimal MC contri- 
bution. All these aspects will be further discussed in 
a forthcoming paper extending the sample to all the 
outliers analysed in Rll and all the local (U)LIRGs of 
Vega et al. (2008). 

A bimodality in the star formation efficiency 
(SFE=SFR/M gas ) and Kennicutt-Schmidt relation 
(Ssfr oc ) of local and high-z mergers and non- 
or weakly-interacting star forming galaxies based on 
CO measurements, has been recently suggested by 
Gcnzel et al. (2010) and Daddi et al. (2010). Our 
upper and lower envelopes in Fig. 3, within which all 
our sources are contained, may correspond to these 
two sequences. Further investigation with larger and 
unbiased samples of galaxies will help in this sense. 
The latter would be explained by the fact that, while 
locally galaxies with IR luminosities exceeding 10 12 
Lq are predominantly associated with merging events 
(e.g. Dasyraetal. 2006), at high-z these high- IR lumi- 
nosities may be achieved mostly via cold gas accretion 
(Powell et al. 2011). Our study confirms that high-z 
(U)LIRGs are more likely up-scaled versions of normal 
galaxies rather than equivalents of local (U)LIRGs in 
terms of the mode of SF (e.g. Symeonidis et al. 2009, 
Rujopakarn et al. 2011). 

With our current assumptions (see §[3]), all our mid- 
IR selected galaxies appear to include massive popula- 
tions of old (> 1 Gyr) stars and, at the same time, to 
host a moderate ongoing activity of SF (SFRio < 100 
M Q /yr, cf. Fig. [3] (right)). The bulk of stars appear to 
have been formed a few Gyr before the observation in 
essentially all cases. Average estimates can be inferred 
from Fig □ for the z~ 2 (U)LIRGs having the SFH 
shown in the central panel of Fig. [IJ about 66-80% of 
the stellar mass has formed after ~ 1 Gyr from the be- 
ginning of their star formation activity, proportionally 
higher at higher SF efficiency. For the z~2 (U)LIRGs 
presenting, instead, more regular SFH, as those shown 
in the bottom panel of Fig. 02 about 30-43% of stellar 
mass is formed after 1 Gyr. Finally for the LIRGs, on 
average, ~ 28% of the stellar mass has already been 
formed after the first 1 Gyr. 



5. Stellar mass determinations 

Another key physical parameter of high-rcdshift 
galaxies is the stellar mass M*. Our best-fit estimates 
appear in Table [2] where we compare our best-fit M* 
to those derived by F10. They fit BC03 models to the 
optical-to-8 fj,m SED with the HYPERZ code (Bol- 
zonella et al. 2000) allowing different SFHs, from a 
constant star formation to an exponentially declining 
SF with different e-folding times (r- model), and as- 
suming the Calzetti et al. (2000) attenuation law to 
account for dust extinction effects. The so called r- 
model is essentially an analytical approximation of the 
Schmidt's law ip{t) oc exp(— t/r), with T=\/vs c h.- As 
shown in Tab. 03 GRASIL M* are systematically larger 
than those based on HYPERZ. This discrepancy is 
higher for the ULIRGs at z ~ 2, for which it reaches 
a factor of ~ 6 in M* (median factor 2.5), while it is 
lower for the z ~ 1 LIRGs whose median factor is 1.4. 

5.1. Origin of the stellar mass discrepancy 

Michalowski et al. (2012, M12 hereafter) found the 
adopted specific SFH (double component versus r- 
model) to be the major factor affecting the derived 
stellar masses of high-z and highly star forming SMGs, 
accounting for a factor of ~ 2.5 in M*. We have in- 
vestigated the possible origin of the larger differences 
in our M* estimates by considering different models 
for the SFH, continuous vs starburst. In both cases 
we found similar M* discrepancies with respect to the 
HYPERZ ones. As discussed in detail in M12, the 
different cvoutionary models may also contribute to a 
factor 2.5 difference in M* estimates when a r- model is 
considered. We tested this hypothesis by modeling our 
(U)LIRGs with similar SFHs (tau-models) , age, metal- 
licity and Ay plus Calzetti attenuation law and by fit- 
ting the optical data alone as in the HYPERZ code. 
By adopting similar prescriptions we found our M 4 to 
be in full agreement, (i.e. within 10%) with the HY- 
PERZ estimates. Taking into account the differences 
in the specific parameters ruling the SFH (e.g. r) and 
on the SSP libraries used (BC03 and GRASIL SSPs 
rely on the same Padova isochrones but the latter in- 
cludes dusty envelopes around AGB stars) we can con- 
clude that the different SFH (i.e. different r) and SSPs 
cannot account for the larger discrepancies, we find 
among the more obscured ULIRGs at z ~ 2; for these 
the dominant factor is the dust exinction. Of course 
the effect of the different SFH on the stellar mass esti- 
mates becomes larger when comparing r-modcls with 
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Fig. 6. — Top-Left: Correlation between the difference in stellar mass estimates based on the GRASIL and HYPERZ 
codes and the average extinction Ay- In the inset the LIRGs values are shown. Top- Right: \ 2 values for spectral 
models of a typical z~ 2 ULIRG are colour-coded (% 2 increasing from dark to light colours) as a function of the 
average extinction and the stellar mass (see Tab. [5]). Here we show all the solutions with x 2 < 100 (gray dots) and 
with x 2 within Xmin+^O (colored points). Bottom: vs M2dust (from Santini et al. (2010)). Blue squares are 
local ULIRGs, red triangles refer to local spirals and black circles correspond to high-z SMGs. Red filled circles and 
cyan filled squares represent, respectively, our z^ 2 and 1 (U)LIRGs whose stellar and dust masses have been 
rescaled to the Scalo IMF adopted by Santini et al. Median ler error bars are shown on the left. Solid cyan and 
dashed (dot-dashed) orange lines are the predictions of Calura et al. (2008) model for spirals and proto-ellipticals 
with mass of 10 11 (10 12 ) M . 
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Fig. 7. — Left:Thc Lboi estimated from the UV to FIR fit by GRASIL is almost a factor of 4-5 higher than the 
Lboi given by the optical only fit, at least for those cases were we found large differences in stellar mass estimates. 
For this galaxy we estimate a M* ~ 4 times higher. Right: Here the effect is smaller. The average extinctions are 
similar and also the bolometric luminosities. For details see tab. [2] 



SFH characterized by recent (last 50 Myr) burst of star 
formation, as shown by M12. All these galaxies are, 
however, normal star forming galaxies and the only five 
objects requiring the burst on top of the Schimdt-type 
SF are all observed as late-SBs. 

5.2. Extinction as the main source of mass dis- 
crepancy 

Using the GRASIL best-fits, we estimate the av- 
erage total extinction in the rest-frame U-band, Ay, 
from the ratio of the extinguished (Ly) to the unex- 
tinguished starlight (Ly): 

A v = -2.5log{L v /L v ) (2) 

The average Ay values, reported in Table [2 are sys- 
tematically higher for the z ~ 2 than for the z ~ 1 ob- 
jects. These differences are traceable to the different 
evolutionary stages at which the two classes of sources 
are observed. Our best-fit Ay are also systematically 
higher than those obtained with HYPERZ using the 
Calzetti approximation in combination with the fore- 
ground screen of dust, expecially for the z ~ 2 ULIRGs 
for which the median Ay is 2.45 and 1.40 for the two 



cases, respectively. 

As shown in fig. [5] (top-left) wc find, for our 
(U)LIRGs, a tight correlation (r ~ 0.86 for z ~ 2 
(U)LIRGs and r ~ 0.96 for z ~ 1 LIRGs) between 
the difference in the M* estimates based on GRASIL 
and HYPERZ solutions and the average extinction 
Ay , with the difference increasing as a function of the 
extinction Ay. This clarifies why for ULIRGs, with 
higher Ay on average, we see the largest stellar mass 
discrepancies. The stellar mass missed by fitting opti- 
cal data alone is hidden in dust. In concurrence with 
our results, Wuyts et al. (2009) also found that the 
mass underestimate is more severe during the dusty, 
peak SF phase. 

It is worth noting that the availability of the full 
wavelength coverage also plays a crucial role in con- 
straining the models. In fact the first basic constraint 
that is provided by the availability of UV-to-sub-mm 
data is the total Lboi of the galaxy, which otherwise 
should be guessed. This is easier when the optical 
contains most of the total energy. The shape of the 
optical- UV SED indicates some level of dust, but there 
may be stellar populations totally obscured by dust 
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(typically the newly born), and even less young stellar 
populations may be partly extinguished. Without a 
constraint on the optical-UV extinction, i.e. on how 
much is re-emitted in the IR, it is difficult to quan- 
tify how much stellar light (therefore mass) is missing 
from the optical. Two examples are given in Fig. [7] 
which compares the total Lboi inferred from the un- 
extinguished starlight component computed with HY- 
PERZ (in red) to that one computed with GRASIL 
(magenta long-dashed line) and to the total Lir in- 
ferred from the GRASIL best-fits. The two objects 
shown in this figure represent two extreme cases, the 
one on the left corresponding to a large discrepancy 
in stellar mass (~ a factor of 4) and the other one 
on the right characterized by the lowest stellar mass 
discrepancy (~ a factor of 1.17). Of course, due to en- 
ergy balance, the energy of the pure stellar SED must 
be equal to that of the processed SED which fit the 
data. Therefore, the IR part of the SED is a further 
fundamental constraint when it contains an important 
fraction of the energy emitted by stellar populations. 

Figure [6] (bottom) shows that our stellar and dust 
masses are consistent with published estimates for 
high-z SMGs at comparable rcdshifts. Compared to 
local spirals and local (U)LIRGs our z ~ 2 (U)LIRGs 
have higher dust-to-stellar mass ratios, by a factor of 
~ 30 and 6, respectively. Similar results were found by 
Santini et al. 2010 for the high-z SMGs. The M du8 t/M* 
of z^ 1 LIRGs are, instead, similar to those of local 
(U)LIRGs. M* and Mdust of galaxies in this work ap- 
pear also to be in agreement with the predictions of 
the Calura et al. (2008) chemical evolution models for 
proto-ellipticals. The distribution of our sources in this 
plot is also similar to that found by Rowan-Robinson 
et al. 2010 for the Hermes Lockman sample. These 
specific aspects will be dealt with in detail in one of 
the forthcoming papers of the series. 

The top-right panel of Fig. [5] illustrates the possi- 
ble degeneracy in the model solutions, between the 
average extinction and stellar mass, color-coded by 
the value of \ 2 , for a typical z ~ 2 ULIRG. As we 
see, among the many solutions considered, acceptable 
best-fits, within Xmin+^O, arc clearly identified in the 
parameter space and not much degeneracy is apparent. 
In particular Ay and M* seem to be well constrained 
within AA V ~ 0.1 and AM* ~ 0.04 dcx. These val- 
ues refer to the case of fixed SFH, that is to say fixed 
Tj n f and z^sch- If we consider in addition to the many 
combinations of GRASIL parameters also the different 
combinations of parameters ruling the SFH, the scat- 



ter around the best-fit solution increases to AAy ~ 
0.2 and AM* ~ 0.15 dex, well within the typical un- 
certainties for this kind of analysis. This is likely to be 
a result of the continuous spectral coverage from UV 
to sub-mm for our sample. Another likely factor is the 
self-consistent radiative transfer modelling, accounting 
for the distribution of dust with respect to young and 
old stars, from which both the age-dependent dust at- 
tenuation and dust reprocessing are constrained. 

6. Summary and Conclusions 

The present work reports on a selection of results 
from a physical self-consistent study of IR luminous 
galaxies at z ~ 1 and 2. For this study we combined 
far-infrared photometry from two related deep surveys, 
PEP and HERMES, carried out with the Herschel 
Space Observatory, with deep photometry in the near- 
IR and mid-IR with the Spitzer IRAC and MIPS im- 
agers, as well as complete mid-infrared spectroscopic 
follow-up with the Spitzer's infrared spectrograph IRS. 

Our working sample includes 31 Luminous and Ul- 
traluminous Infrared Galaxies selected in GOODS-S 
with 24 (im fluxes between 0.2 and 0.5 mJy. Among 
them 10 are at z~l and 21 are at z~2. They make the 
faintest high-rcdshift galaxies observed spectroscopi- 
cally with IRS, with an excellent suite of photometric 
data in all bands from UV to the millimeter, including 
our key far-IR and sub-mm photometry from Herschel. 

These data have prompted us to investigate the na- 
ture and the main physical properties - like stellar 
mass, bolometric luminosity, star-formation history, 
extinction, as well as the mass assembly history - for 
every object. The novelty of our approach consists in 
modeling the data with a self-consistent physical code 
(GRASIL), exploiting a state-of-the-art treatment of 
dust extinction and reprocessing. 

We find that all of our galaxies require massive pop- 
ulations of old (> 1 Gyr) stars and, at the same time, 
host a moderate ongoing activity of SF (SFRio < 100 
M /yr with SFRio being the SFR averaged over the 
last 10 Myr). Our detailed analysis appears to give 
important hints also on the past history of SF with 
the bulk of stars consistent with having been formed a 
few Gyr before the observation in essentially all cases. 
Only five galaxies of the sample seem to require, on the 
basis of an application of the F-test to our x 2 analysis, 
a recent starburst superimposed on a low-level, secu- 
larly evolving star formation history. Even in these 
five objects, however, the gas mass involved in the SB 
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amounts to only a small fraction of the galactic mass, 
< 4%, and all of them are observed just at the end of 
the burst event (i.e. late SBs). 

We find substantial discrepancies between our re- 
sults and those based on optical-only SED fitting for 
the same objects. Our physically consistent best-fit 
model solutions of the observed SEDs indicate higher 
extinctions (by AAy ~ 0.81 and 1.14) and higher stel- 
lar masses (by ALog(M^) ~ 0.16 and 0.36 dex) for 
z~l and z^2 (U)LIRGs, respectively. The stellar mass 
difference is larger for the most obscured objects and 
correlates with the total dust extinction. This is in 
agreement with the results of Wuyts et al. (2009) who 
also found that the mass underestimate is more severe 
during the dusty, peak SF phase. 

We also find lower SFRs on average than those com- 
puted from the total Lm, using the Kennicutt (1998) 
calibration due to the significant contribution to the 
dust heating by intermediate-age stellar populations 
through 'cirrus' emission: only ^27% and ^34% of 
the total Ltr appears to be due to ongoing SF in MCs 
for z ~ 1 and z ~ 2 (U)LIRGs, respectively. 

We have demonstrated with the present study the 
importance of a self-consistent treatment of dust ex- 
tinction and reprocessing effects in luminous star- 
forming galaxies via radiative transfer modelling. In 
future works many aspects will be further investigated, 
such as a more complete understanding of the degen- 
eracies in our model solutions, including exploring the 
effects of also other IMFs. This is of particular rel- 
evance for a further refinement of the galaxy star- 
formation histories and the estimates of stellar mass 
functions and SFR functions, that have a significant 
impact on the studies of galaxy formation and evolu- 
tion. 
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APPENDIX: Best-fit SEDs of the entire sample of high-z (U)LIRGs 
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Fig. 8.— GRASIL best-fits (solid black line) to the observed SEDs (red circles) of z ~ 1 and z ~ 2 (U)LIRGs. 
IRS spectra appear in the inset window (magenta line). The color-coded lines represent the unextinguished starlight 
(orange dot-dashed), extinguished starlight (green dashed), cirrus emission (magenta dotted) and MC emission (blue 
long dashed). 
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Fig. 9. — same as in fig. [S] 
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Fig. 10. — same as in fig. [5] 
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